%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% COMPUTE AVOIDED EMISSIONS OF CO2

% "DEFORESTATION IN THE AMAZON:
% A UNIFIED FRAMEWORK FOR ESTIMATION AND POLICY ANALYSIS"

% by Eduardo Souza-Rodrigues

% This version: November 2018

% OBSERVATIONS: 
% 1. Run this program after running "dem_def_reg.m"
% 2. The data set has to be loaded before running this program.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% (A) COMPUTE SUPPLY OF AVOIDED CARBON EMISSIONS

% I. Compute supply only for the main specification
if prod_index == 2  % Crops-Pasture Index 
    
% II. Prepare the Data

% Set vector of carbon prices, US$/tCO2
co2_price1 = 0.1:0.1:0.9;
co2_price2 = 1:0.5:35;
co2_price  = [0 co2_price1 co2_price2];
clear co2_price1 co2_price2

% NOTE: Convert Difference of Carbon Stock between Forest and Deforested to CO2 Emissions
% CarbonDiff       = CarbonStock_Forest - CarbonStock_Deforested (tC/ha)
% CO2_Emissions    = (44/12) * CarbonDiff,
% because:    X tC = (44/12) * Y tCO2
% So, I use:
% co2_emissions_for = (44/12) * CarbonDiff; 


% III. Compute Counterfactual Emissions of Carbon

% Store Results
forest_price_ha   = -10*ones(T,size(co2_price,2));  % Local price US$/ha considering difference of carbon stock in Forest and Deforested areas
SI_ivqr           = zeros(T,size(co2_price,2));
S_ivqr            = zeros(T,size(co2_price,2));
Deforest_ivqr     = zeros(T,size(co2_price,2));
avoid_carbon_ivqr = zeros(T,size(co2_price,2));
tot_carbon_ivqr   = zeros(size(co2_price));


% Fix Carbon Price
for p = 1:size(co2_price,2)
    
    % Fix Observation (= Municipality)
    for j = 1:T
                       
        % 1. Convert prices from US$/tCO2 to US$/ha
        forest_price_ha(j,p) = co2_price(p) * (44/12) * CarbonDiff(j);
        % NOTE: US$/ha prices are local because of differences in carbon stock across municipalities
        
        % 2. Counterfactual Deforestation for Price = US$/ha (where US$ 1 = R$ 2.17, in 2006)
        if beta_ivqr(1,Uivqr(j,1)) < 0
            
            % Single-Index
            if penalty == 0     % Actual Monitoring (Main Specification)
                
                SI_ivqr(j,p)  = Xc(j,:) * beta_ivqr(:,Uivqr(j,1)) + beta_ivqr(1,Uivqr(j,1)) * (forest_price_ha(j,p)/Qa(j,1));
                
            elseif penalty == 1 % Fines post-2004 at pre-2004 level (i.e., Fines 2004-05 = Fines 2003-04)
                
                SI_ivqr(j,p)  = Xc(j,:) * beta_ivqr(:,Uivqr(j,1)) + beta_ivqr(1,Uivqr(j,1)) * (forest_price_ha(j,p)/Qa(j,1))...
                    - beta_ivqr(2,Uivqr(j,1)) * Fines2005(j) + beta_ivqr(2,Uivqr(j,1)) * 2 * Fines2003(j);
                
            end
            
            % Share of Land-use
            S_ivqr(j,p) = (exp(SI_ivqr(j,p)) / (1 + exp(SI_ivqr(j,p))));
            
            % Deforested
            Deforest_ivqr(j,p) = S_ivqr(j,p) * Area(j,1);
            
        end
        
        % 3. Counterfactual Avoided Emissions of Carbon        
        avoid_carbon_ivqr(j,p) = CarbonDiff(j)    * (Deforest_ivqr(j,1) - Deforest_ivqr(j,p));
        % NOTE: Avoided Emissions  = Emissions per ha * Deforested Area in hectares
                  
    end
    
    % 4. Total Counterfactual Avoided Emissions
    tot_carbon_ivqr(p) = sum(avoid_carbon_ivqr(:,p) .*(CarbonDiff > 0) .*(beta_ivqr(1,Uivqr)<=0)');
               
end



%% (B) SAVE RESULTS

if farm == 1
    
    if penalty == 0         
        
        % Total Avoided Emissions
        Total_emissions_farm_1 = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_1   Total_emissions_farm_1   -ascii
        
    elseif penalty == 1
        
        % Total Avoided Emissions
        Total_emissions_farm_1_fines = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_1_fines   Total_emissions_farm_1_fines   -ascii                
        
    end        
    
elseif farm == 2
    
    if penalty == 0         
        
        % Total Avoided Emissions
        Total_emissions_farm_2 = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_2   Total_emissions_farm_2   -ascii
        
    elseif penalty == 1
        
        % Total Avoided Emissions
        Total_emissions_farm_2_fines = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_2_fines   Total_emissions_farm_2_fines   -ascii
        
    end    
        
elseif farm == 3
    
    if penalty == 0         
        
        % Total Avoided Emissions
        Total_emissions_farm_3 = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_3   Total_emissions_farm_3   -ascii
        
    elseif penalty == 1
        
        % Total Avoided Emissions
        Total_emissions_farm_3_fines = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_3_fines   Total_emissions_farm_3_fines   -ascii
        
    end       
   
elseif farm == 4
    
    if penalty == 0         
        
        % Total Avoided Emissions
        Total_emissions_farm_4 = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_4   Total_emissions_farm_4   -ascii
        
    elseif penalty == 1
        
        % Total Avoided Emissions
        Total_emissions_farm_4_fines = [co2_price' tot_carbon_ivqr'];

        % Save
        save Total_emissions_farm_4_fines   Total_emissions_farm_4_fines   -ascii
        
    end
        
end

end

%% (C) DELETE IF NOT NEEDED

clear p j
clear forest_price_ha SI_ivqr S_ivqr Deforest_ivqr
clear avoid_carbon_ivqr tot_carbon_ivqr


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
